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ABSTRACT 

Observations indicate that the gaseous circumstenar disks around young stars vary significantly in 
size, ranging from tens to thousands of AU. Models of planet formation depend critically upon the 
properties of these primordial disks, yet in general it is impossible to connect an existing planetary 
system with an observed disk. We present a method by which we can constrain the size of our own 
protosolar nebula using the properties of the small body reservoirs in the solar system. In standard 
planet formation theory, after Jupiter and Saturn formed they scattered a significant number of 
remnant planetesimals into highly eccentric orbits. In this paper, we show that if there had been a 
massive, extended protoplanetary disk at that time, then the disk would have excited Kozai oscillations 
in some of the scattered objects, driving them into high-inclination {i > 50°), low-eccentricity orbits 
{q > 30 AU). The dissipation of the gaseous disk would strand a subset of objects in these high- 
inclination orbits; orbits that are stable on Gyr time scales. To date, surveys have not detected 
any Kuiper Belt Objects with orbits consistent with this dynamical mechanism. Using these non- 
detections by the Deep Ecliptic Survey (DES) and the Palomar Distant Solar System Survey we are 
able to rule out an extended gaseous protoplanetary disk {R^ ^80 AU) in our solar system at the 
time of Jupiter's formation. Future deep all sky surveys such as the Large Synoptic Survey Telescope 
(LSST) will all us to further constrain the size of the protoplanetary disk. 



1. INTRODUCTION 

Understanding the properties of the gas rich disks 
around young stars is crucial to understanding the for- 
mation of planetary systems. Observations demonstrate 
that these protoplanetary disks exhibit sig nificant diver- 
sity in their sizes, masses, and lifeti mes fe.g. Ha isch et al.l 
HOOI; Mamajek 2009; An drews et a l. 2009, 20lJ). These 
variations must affect the ubiquity and diversity of ob- 
served exoplanet systems. However, it impossible to con- 
nect any of the Gyr-old existing planetary systems di- 
rectly to the observed few Myr-old disks. However, our 
own solar system provides a unique test case. In this one 
system we have access to a significant amount of data, 
stored in the small body reservoirs, that allows us to 
make inferences about the early state of our own solar 
nebula. 

There is a long tradition of using observations from 
the present day solar system to try to constrain the 
properties of the original solar nebula. Early on, re- 
searchers noted that they could define a minimum mass 
solar nebula (MMSN), the lowest mass disk that could 
have created the planets in their current locations, by 
augmenting the observed solids in the planets with hy- 
drog^en and helium unt i l they reach s olar composition 
(|Weidenschillind Il977bl : iHavashil IT98lh . The resulting 
power-law profile has served the standard benchmark 
model in planet formation simulations. More recently, 
researchers have suggested that an a pparent hard edge 
to the classical Kuiper belt at 50 AU (| Allen et al.ll2001f ) 
may be indicative of the initial size of the protoplanetary 
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disk. Indeed, lAdamsl (|2010f ) and references within have 
used this apparent outer edge to the solar system as part 
of an effort to understand the birth environment of our 
solar system. For example, Adams presents arguments 
that a star cluster dense enough to truncate the disk at 
around 50 AU is consistent with a cluster capable of pro- 
ducing a super-nova able to contaminate our early solar 
system with short lived radioactive isotopes. Although 
this paints a self-consistent picture, it is important to 
note that the mechanisms that truncate the outer re- 
gions of the disk (especially photoevaporation) actually 
act primarily on the disk gas^ while the data (the dis- 
tribution of planets and KBOs) refers to the disk solids. 
There are well-known physical me chanisms capable of de- 
coup ling the solids a nd the gas (jWeidenschillind Il977al : 
i.Youdin Shul[2002h . Indeed, some authors have even 
argued that a large reservoir of drifting solids born in an 
extended gaseous disk could be useful in enhancing the 
solid to gas ratio, increasing the likeli hood of planetesi- 
mal f ormation and core formation (e.g. lYoudin fc ChiangI 
[200l . Therefore one cannot claim that an "edge" in the 
planetesimal disk necessarily corresponds to an edge to 
the gas disk. 

In this paper we present an argument that constrains 
the size of the gaseous disk at the time of Jupiter's forma- 
tion. Planet formation is an inefficient process; current 
models suggest that the total mass of planetesimals left 
over after the planets are fully assembled is equivalent 
to, or up to 5 times greater than, the total mass of solids 
incorporated into the plane ts (e.g. "Pollack et al.'"1996'; 
Ussauer 1987; Tho mmes et al. 2003; Alibert et al. 2005' 
iHubickvi et al.l 120051 and see iThommes Duncanl 120061 
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for a review) . Therefore during and shortly after the for- 
mation of the giant planets, must have been a significant 
number of planetesimals left over in the planets' feeding 
zones. These planetesimals are fated to undergo close en- 
counters with the planets, and many will be scattered to 
high eccentricity orbits. In a simple 3-body scattering, 
the perihelion of the scattered object remains roughly 
unchanged, meaning that the planetesimal will return to 
undergo additional close encounters. Eventually these 
objects will either collide with the planet or be ejected 
by the system. 

However, in order to form the planets there must also 
have been a circumstellar gas disk, a disk that poten- 
tially could have extended beyond the orbits of the plan- 
ets. And the presence of a massive gaseous disk can 
alter these scattered planetesimals' orbits. For planetes- 
imals between a few tens to a few thousands of kilometers 
in si ze, they are too lar ge to experience significant gas 
drag (|Adachi et al]|1976f ), and yet are too small to un- 
dergo any sort of tidal interaction with the disk (|Wardl 
Il997l and references within). Given the cylindrical sym- 
metry of the disk, the major gravitational effect of the 
disk is of the Kozai type (Kozai 1962), which arises when 
a perturber can be treated as a circular ring. In the 
presence of the potential from a massive gaseous disk, 
planetesimals scattered from the Jupiter/Saturn region 
into moderate inclination, high eccentricity orbits may 
undergo Kozai oscillations. 

During a Kozai oscillation, the presence of perturbing 
material exterior to the scattered particle forces the par- 
ticle exchange eccentricity for inclination, all the while 
conserving the z-component of the angular momentum. 
In the situation presented here the gaseous disk serves 
as the perturbing potential. The scattered object os- 
cillates between a low-inclination, eccentric orbit and an 
high-inclination, more circular orbit. During the high ec- 
centricity phase the particle is susceptible to additional 
close encounters with the planets, but during the low ec- 
centricity phase, the particle's orbit will no longer cross 
the planets. After the disk disperses, statistically some 
of the objects will be stranded at high inclination orbits. 
This dynamical mechanism is capable of producing ob- 
jects with large enough perihelia to place them in Kuiper 
Belt on orbits that are stable for the lifetime of the so- 
lar system with the planets in their present day orbital 
configuration. But while these particles have semi-major 
axes that would place them in the Kuiper Belt, they 
would have substantially higher inclinations than what 
has been observed so far. In this paper we explore this 
dynamical mechanism as a way to use the current dis- 
tribution of Kuiper Belt Objects to place constraints on 
the early structure of the protoplanetary disk. 

In 0we describe our numerical method and model pa- 
rameters. In ^we discuss the outcomes for various disk 
parameters. In g]we discuss how we can use these re- 
sults to constrain the size of the protosolar nebula using 
resul t s from the Deep Ecliptic Survey (DE S ; iMillis et al.l 
[2QQllElliot et al.ll2QQ5l : iGulbis et a l. 2010) and the Palo- 
mar Distant Solar System Survey (Schwa mb et al.ll2QQ9L 
[2010). And finally, in g5] we summarize our results and 
discuss the implications. 

2. MODEL 



In this paper we investigate the behavior of planetes- 
imals initially in the giant planet regions and model 
their gravitational interactions with the planets and the 
gaseous disk. We initialize these simulations with four 
giant planets in the compact configuration proposed by 
iTsiganis et al. (2005,) as the initial positions of the plan- 
ets in the Nice Models with Jupiter, Saturn, Uranus and 
Neptune at 5.45, 8.18, 11.5, and 14.2 AU, respectively. 
However, as we discuss in ^ the results of this pa- 
per are insensitive to this choice. For each simulation 
we calculate the orbital evolution of TVtp = 6 x 10^- 
10^ massless test particles initially distributed from 4- 
15 AU with small eccentricities (e = 0.1) and incli- 
nations (0° < z < 3°) to represent the planetesimals 
left over in the planet forming region. We calculate 
the orbital evolution of the scattered planetesimals us- 
ing the SWIFT RMVS3 numerical integrator package 
(jLevison fc DuncanI Il994l ). We integrate the orbits for 
5 Myr with a time- step of 0.15 yr, and we remove par- 
ticles from the simulation if they become unbound or 
they experiences a collision with a planet or the Sun. 
Additionally, as we are interested only in particles whose 
orbits are not perturbed by the potential due to the birth 
cluster or the galaxy, we remove particles if their semi- 
major axis (a) becomes larger than 300 AU. 

To calculate the particles' gravitational interaction 
with the disk gas we numerically integrate Poisson's 
equation over a cylindrically symmetric grid. This allows 
us to obtain the potential and the accelerations derived 
from the disk gas on the particles using bicubic inter- 
polation. For computational simplicity we use simple a 
simple parametric form for the disk gas. We describe the 
surface density profile as a function of distance from the 
sun r at time t as a truncated power-law. 



S(r,t) 



6(r, i?D)exp ( — 



(1) 



where Q{r^RD) describes the functional form of the trun- 
cation of the outer disk at size Rd^ and td is the disk 
depletion timescale. For the fiducial disk models we 
use Eo = 2000g cm~^ and 7 = 3/2 cons istent with th e 
benchmark minimum mass nebula profile fflavas hilll981f ) . 
We investigate the effect of the sharpness of the disk trun- 
cation on distribution of planetesimals by comparing two 
truncation profiles, a sharp cutoff. 



e{r,RD) = ecut{r,RD) 



1, if 1 AU < r < 
0, otherwise. 



and an exponential cutoff 

exp(-r/i?i)), if r > 1 AU, 



0exp(r, Rd) 



0, otherwise. 



(2) 



(3) 



(c.f. [Andrews et"aIll2009L I2010D . In ah disk models we 
truncate the inner disk at 1 AU, but we vary the outer 
disk cutoff {Rd) from 30 AU to 100 AU. To preserve the 
stability of the integration we use a "puffy" exponential 
profile for the disk vertical profile. 



p{r,z) = E(r)exp (-^) 



(4) 



where H is constant throughout the disk. This disk 
model is chosen for its simplicity and is not intended 
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to represent a detailed disk evolution model. 

There are five free parameters for the disks in these 
simulations: Sq, ^d^ Rd^ and 7, and the functional form 
of 6 (r,i?D). In total we present the results of the evo- 
lution in 35 different gas disks. Our fiducial parameters 
are Eq = 2000 g cm-^, = 2 Myr, 6(r,i?i)) = Ocut, 
Rd = 100 AU, and 7 = 3/2. With the fiducial disk pa- 
rameters, we perform 7 simulations varying from 20 
to 4000 g cm~^. We also perform 6 simulations varying 
td from 2 x 10^ to 4 x 10^ years. Using the fiducial Eq 
and Tjj we then perform simulations with Rjj =30, 40, 
50, 60, 80 and 100 AU. For the sharp cut off we allow 
7 = 3/2 (the MMSN value) and 7 = 1 (a value consistent 
with some viscous evolution models). For the exponen- 
tial cutoff we us e 7 = 1/2 and 1, app ropriate values for 
observed disks ([Andrews et al.|[2010f ). In addition, we 
compute one simulation in the fiducial disk with the four 
giant planets in the configuration of the modern day solar 
system. We run these scattering simulations for 5 Myr 
(10 Myr for the disk with = 4 x 10^). 

In most of these simulations, some particles are cap- 
tured by the Kozai resonance and forced into orbits with 
q > 30 AU; orbits that are potentially stable in the mod- 
ern day solar system. However, although a large number 
of particles undergo oscillations with large perihelion val- 
ues at some point during their orbital evolution, at the 
end of the scattering simulations only a handful of these 
particles happen to have been caught on the q phase of 
their orbits. Still, the time at which a given particle is 
scattered into one of these orbit is random, and could 
as easily have happened right before the disk dispersed 
as at any other time. Therefore instead of restricting 
our stability calculation to only those particles still ac- 
tive at the end of the simulations, we randomly extract 
the orbital elements (a,e, and i) of particles with large 
periheha {q > 20 AU) values throughout the simulations. 
We then randomly assign orbital angles and create a new 
suite of equally probable, potentially stable particles. We 
integrate these new particles in the presence of the four 
giant planets, in their modern day configuration, for 4.5 
Gyr. We use these particles to define the probability 
that a particle will survive as a function of its perihe- 
lion distance at the beginning of the stability calculations 
{Ps{q)). We can then return to the initial scattering cal- 
culations and use the probability that each particle (at 
its given perihelion distance) will be stable to calculate 
the expected fractional retention of particles, 

f = J2 Ps{qi)/NTP. (5) 

1=1 

Here Ps{q) is the function defined by the stability cal- 
culation, and qi refer to the perihelion distance of each 
particle at the end of the initial scattering simulations. 
Because Ps{q) goes to zero as qi approaches 20 AU in 
all simulations, we are confident that the initial cut-off 
of g = 20 AU does not impact the statistics in any part 
of this study. To determine our 95% confidence inter- 
vals, we assume that our data is distributed according 
to a binomial distribution and u se a non-inform ative, 
Bayes-Laplace uniform prior (see iCameroril 1201 IL for a 
nice description of this method). 

In these simulations, the particles only interact with 



the gas though its gravitational p otential, not through 
aerodynamic drag. As indicated in iBrasser et al.l (|2007l ) 
the gas drag is only dynamically important for objects 
with diameters smaller than about D = 30 km. The im- 
pact of the gas drag on these small bodies is to circularize 
the orbits, and thus they do not participate in the effect 
explored in this paper. Furthermore, as we discuss in §31 
most of the mass is likely to be in 100 km, or larger, sized 
objects and smaller objects would not be detectable in 
current datasets. Therefore we only concern ourselves 
with these larger bodies and neglect gas drag in these 
simulations. 

3. RESULTS 

In figure [1] we compare the orbital distribution of plan- 
etesimals produced by disks of various sizes after the scat- 
tering simulations, by which time the gas disks have dis- 
persed. Without the gravitational potential of a gaseous 
disk (left-hand panels), the planetesimals are scattered 
by the planets into eccentric and moderately inclined or- 
bits, but the perihelia (g) remain roughly unchanged. 
Therefore, all of the particles remain on planet-crossing 
orbits that are unstable on timescales short compared to 
the age of the solar system. These unstable particles are 
indicated by gray points. In contrast, in middle and right 
panels we show the particles' evolution in the presence of 
a disk with a MMSN profile (7 = 3/2) and a sharp cutoff 
(9 = 6cut), truncated at Rd = 30 AU and 100 AU, re- 
spectively. In both of these simulations, a subset of the 
particles are placed into orbits with larger perihelia and 
higher inclinations. In the small disk, some particles have 
their perihelia raised so that they have a high probability 
of being stable in the early compact solar system, with 
q > 20 AU (particles marked in blue), but it is only in 
the larger disk that the perihelia are raised to the point 
where they would be stable in the modern solar system, 
q > 35 AU (marked in red). 

These particles get onto these high inclination, low ec- 
centricity orbits via Kozai oscillations. In figure [2] we 
show the temporal evolution of an example particle in 
the 100 AU disk. The particle undergoes a series of close 
encounters with Jupiter and is scattered out to a high ec- 
centricity orbit. After about half a million years we see 
a number of intervals in which the perihelion distance 
increases and decreases in phase with oscillation in the 
inclination while the semi-major axis remains constant, 
as is characteristic for Kozai oscillations. As further ev- 
idence, in the right hand panel we plot the eccentricity 
and the argument of perihelion for this particle. The 
large circles indicate a time interval (shaded in the left 
hand panel) of a single Kozai cycle. In this time the argu- 
ment of perihelion oscillates, indicating that this is a true 
Kozai cycle. Additionally, by comparing figure [2t with 
the Kozai dynamics exerted by t he planets on objects 
with larger semi- major axes (c.f. iThomas fc Morbidellil 
[19961 ) we can be confident that the perturbations seen 
here are dominated by the external disk potential not 
perturbations from the interior planets. During each of 
the periods of high eccentricity (small q) the particle un- 
dergoes additional encounters with Jupiter and Saturn, 
causing sharp jumps in the particle's semi-major axis. As 
the disk is depleted the timescale for the oscillation be- 
come longer. In the end the period of oscillation becomes 
infinite, stranding the particle on a high inclination orbit. 
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Fig. 1. — The eccentricity (top), inclinations (middle), and per- 
ihelion distance (bottom) for particles after 5 Myr in simulations 
with no disk potential (left), and the potential of a MMSN disk 
truncated at 30 AU disk (middle) and 100 AU (right). In blue we 
show particles stable in the compact early solar system and in dark 
red we show particles stable in the present day solar system, while 
test particles on unstable orbits are in gray. 

Regardless of the exact initial orbits of the planets, 
at some point they must have moved into their cur- 
rently observed configuration. Therefore we are really 
only concerned with the orbital distribution of particles 
that would be observable today due the fact that their 
orbits are stable in the present day solar system, for the 
lifetime of the solar system. The location of all stable ob- 
jects produced in simulations with various disk profiles 
and sizes are shown in figure [3l We show four panels, a 
"shallow" and "steep" surface density profile for each of 
our two truncation functions. For the upper plots (sharp 
cutoff) the disk truncation function is a step function 
(B = Ocut)- The surface density profiles are, for panel 
a, the benchmark MMSN 7 = 3/2 and, for panel b) the 
shallower 7 = 1 more consistent with viscous evolution. 
For the lower plots (exponential cutoff) the disk trunca- 
tion function is O = Oexp- We choose to show shallower 
surface density indices, 7 = 1/2 and 7 = 1 in panels c 
and d, respectively, to be co nsistent with observations of 
disks () Andrews et al.ll2QlQf ). The points are color-coded 
by the disk size, with Rp = 40, 50, 60, 80, and 100 AU 
in magenta stars, red diamonds, yellow squares, green 
hexagons, and blue circles, respectively. 

The orbital distribution of stable test particles are 
qualitatively different depending on the size and cut- 
off profile of the disk. In disks with sharp cutoff 
there is a limited semi- major axis range in which stable 
high-inclination particles can be produced. An exterior 
quadrapole moment is much more effective in raising the 
Kozai resonance than an inner on e, and the dynamics of 
the t wo are significantly different (jThomas fc Morbidellil 
|1996[ ). Therefore larger disks can excite Kozai oscilla- 
tions in particles with larger semi-major axes. Addition- 
ally, these disks are well approximated by an annulus 
at the outer edge as that is where the mass is concen- 
trated. The shallower surface density profile (7 = 1) is 
slightly better at producing high inclination particles as 
more of the disk mass is at this outer radii. In disks 
with exponential profiles there is mass at a large range 
of radii. This allows particles to undergo Kozai oscilla- 



tions at large a, however the fact that the disk is not 
dominated by a single annulus of material at the outer 
edge means that the oscillations are shallower. Therefore 
the objects produced do not have very low eccentricity 
orbits. The shallower surface density profile (7 = 1/2) 
again more effectively generates these particles due to 
the fact that there is more mass at larger radii to excite 
Kozai oscillations. 

The effect of disk size on retention can be seen more 
clearly in figure H] in which we show the fractional reten- 
tion (/, defined in eq. [5]) as a function of the disk size 
and profile (using the parametric form shown in equation 
([1])). The downwards arrows indicate the upper limit on 
simulations in which no particle ever has its perihelia 
raised to a potentially stable orbit. There is a trend that 
larger disks tend to retain a larger fraction of particles. 

We note that in these various simulations we main- 
tained the same normalization of the surface density pro- 
file of the gas disk at 5 AU = 180g cm"^). How- 
ever, this arbitrary choice does not significantly impact 
the results. In figure [5] we show the fraction of parti- 
cles retained as a function of disk mass (with the fiducial 
depletion timescale of 2 Myr) and depletion timescale 
(with the fiducial disk mass of 0.05 Mq). Other than at 
very small disk masses and very short disk lifetimes the 
fraction of retained particles in remarkable insensitive to 
these disk parameters. In order to understand this insen- 
sitivity, we examine the timescale for a Kozai oscillation. 
For a « the characteristic timescale for Kozai oscil- 
lations due to a perturbation of a narrow, infinitesimally 
thin annulus of radius R and mass M on a particle or 
semi-major axis a is 



2 Mr 







Stt M 



a J 



3/2 



(6) 



(e.g. iKiseleva et al.l (|1998f )). Despite the fact that the 
systems considered in this paper are more aptly described 
as a series of concentric rings with finite vertical thick- 
ness, the Kozai timescale for the disks with a O = Ocut 
can be reasonably well approximated by modeling the 
disk as an thin annulus of mass M = Mjj and radius 
R = Rd- This timescale becomes longer diS Rd increases, 
and as Mjj decreases. Therefore there is a limiting mass 
for which a single oscillation can occur within the disk 
lifetime. For disks with a timescale for dispersal of less 
than this we do not see any particles excited into orbits of 
interest. As mentioned previously, for our fiducial disk 
of a MMSN truncated at 100 AU, this estimate of the 
Kozai timescale yields 2 x 10^ years. Additionally, for 
the fiducial disk lifetime of 2 Myr, the minimum mass 
required in order to allow for a single Kozai oscillation is 
2 X 10~^ Mq. This means that for our fiducial disk pa- 
rameters, the particles undergo many Kozai cycles, and 
during each close encounter with a planet the particle 
may either be scattered into or out of the Kozai cycle. 
As the disk mass or lifetime increases, a larger fraction 
of the particles undergo Kozai oscillations at some point 
in their evolution. However, most of these particles at 
some point undergo a close encounter with a giant planet 
and are ejected from the system. Therefore the number 
of particles left after disk dispersal asymptotes with in- 
creasing disk mass. So long as Jupiter and Saturn did 
not form so late in the disk evolution that, within a few 
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Fig. 2. — In panel(a) we show the semi-major axis (black curve) and pericenter distance (Grey curve) for an example particle in the 
fiducial simulation. This particle is scattered by Jupiter and then undergoes a series of Kozai oscillations and close encounters. Eventually 
it is stranded in a high inclination orbit. In panel (b) we show the inclination of the particle. The shaded region indicates an example 
Kozai cycle. In panel (c) we show the eccentricity and argument of pericenter in a polar plot for the example particle. The small dots show 
these parameters at all times, while the large circles indicate quantities while the particle is in the shaded region of the left panels. 



tens of thousands of years, the disk mass decreased to 
less than a Jupiter Mass, we are insensitive to the actual 
disk mass and lifetime. This means that we can safely 
assume we are in a parameter space in which disk size 
and profile are important parameters, but disk mass and 
lifetime are not. 

One may also be concerned about the sensitivity of 
these results on our assumption about the position of 
the planets. There are two possible reasons the planetary 
positions may be important, first the secular perturba- 
tions from the interior planets may stabilize the parti- 
cles against the Kozai resonance. Additionally, the loca- 
tion of the planets may scatter the particles in different 
ways. Therefore, we perform the same calculation with 
the planets in the present-day solar system, as opposed 
to the Nice initial conditions. The fractional trapping is 
virtually identical in this case (see the diamond in fig[5|). 
We note that in the fiducial disk model, the total mass of 
the disk gas interior to 30 AU is over 13 Mj^p, therefore 
the secular influence of the planets is minor compared to 
that of the disk. This is because, in order to undergo a 
deep Kozai oscillation, a particle at a given semi-major 
axis must either have a minimum eccentricity or inclina- 
tion. As we are interested in using the Kozai resonance 
to produce high perihelion (low eccentricity) objects, we 
are interested in the initially eccentric objects. In or- 
der for the initially dynamically cold particles need to 
reach these minimum eccentricities, they must undergo 
a close-encounter with a planet. After the scattering the 
particles orbit will have a perihelia close to the orbit of 
the scattering planet {q = a(l — e) ~ a^). Therefore, if 
we look at a given a, a particle scattered by a planet with 
a smaller ap will have a higher eccentricity than a par- 
ticle scattered by a more distant planet. Starting with a 
low-inclination population, the planets beyond ~ lOA/7 
cannot efficiently produce particles with sufficiently high 
eccentricities to undergo Kozai oscillations at the semi- 
major axes of interest. Therefore as long as Jupiter (and 
Saturn, to a lesser extent) were not far from their present 
day locations near the end of the gaseous disk's lifetime. 



the results are not sensitive to the other planets loca- 
tions. In these direct N-body simulation the potentially 
stabilizing impact of the Additionally it is important to 
note that the time t = in these calculations corresponds 
to the time of the formation of the giant planets, not the 
initial time of the formation of the disk. 

4. COMPARISON TO OBSERVATIONS 

In order to convert these retention fractions into a 
form suitable for comparisons with observations, we must 
make a series of assumptions about the properties of 
the initial population of planetesimals, specifically their 
initial abundance, size distribution, and physical char- 
acteristics. Ideally we would like these assumptions to 
be guided by the observed characteristics of populations 
with the same source, unfortunately, we discuss later in 
this section, there is no observed population that we can 
be sure came from the same source region. Therefore 
before we make a detailed comparison to observations 
we will first present an order of magnitude estimate for 
the number of high-inclination particles produced by this 
mechanism. 

Theoretical models suggest that in order to form the 
cores of the giant planets within the observed lifetime of 
the gaseous protoplanetary disks, the initial disk mass 
must have been 2 to 5 or even 10 ti mes larger than the 
minimum mass solar nebula (e.g. Pollack et al. "1996'; 
Lissauer 1987; Thommes et al. 2003; Alibert et al. 2005' 
iHubickvi et all 120051 and see IThommes DuncanI 120061 
for a review) . This means that fewer than half of the ini- 
tial planetesimals were incorporated into the giant plan- 
ets, so they must have been ejected from the solar system 
or scattered into the sun. Therefore, we assume that the 
total amount of mass in our solar system's initial plan- 
etesimal swarm must be at least equal to the mass of 
solids incorporated into the planets, ~ 4OM0, and is 
probably significantly larger. 

It is less clear how this mass was distributed among 
the planetesimals. Without direct evidence for the true 
size distribution of these particles, we assume that the 
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Fig. 3. — The eccentricity and inclination of objects left in disk with disk size = 40, 50, 60, 80, and 100 AU (magenta stars, red 
diamonds, yellow squares, green hexagons, and blue circles, respectively). In the upper panels we shade the regions according to perihelia 
distance (q = 20, 30, 35 AU darker to lighter). 
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Fig. 4. — The fraction of objects left in stable high-inclination 
orbits after 5 Myr as a function disk size and profile for a disk 
with a sharp cutoff (left), and an exponential cutoff (right). The 
circles, squares, and diamonds correspond to 7 = 3/2, 1, and 1/2, 
respectively. The downwards arrows indicate simulations in which 
no particle ever has its perihelia raised to a potentially stable orbit. 



population has a size distribution similar to what is ob- 
served in other small body reservoirs in the outer so- 
lar system-an assumption we will revisit at the end of 
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Fig. 5. — The fraction of objects left in stable high-inclination 
orbits as a function disk mass (left), and disk lifetime (right) for 
a disk with Rd = 100 AU and 7 = 3/2. In the left-hand plot the 
light gray diamond indicates the results for a simulation in which 
the scattering is done by planets on their present day orbits. 

this section. In general it appears that there is a rel- 
atively steep size distribution for large particles, which 
then, at some size 5 < 5: 150 km , rolls-over to a shal- 
lower distribution ([Bernstein et al.l i2004 : .Fuentes et al.l 
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Fig. 6. — The perihelion distance v. the inchnation for all of the 
simulations (blue, green, yellow, and red symbols are the same as 
in figure [3]) compared to the Kuiper Belt objects (+ symbols) and 
centaur/scattered disk objects (x symbols) from the minor planet 
center. 



aser fc Kavelaarsll2QQ9f ) . Particles with sizes near 
this break dominate the mass of the whole population. 
Therefore if we assume a similar size distribution for 
our initial population, this corresponds to greater than 
10^ initial planetesimals of around D = 100 km in size. 
With a capture efficiency around 10~^, we expect the 
mass in our hypothetical reservoir of objects with incli- 
nations larger than 50° to be > O.OSMq, or > 10^ ob- 
jects. This is of the same order as the upper limits on 
the mass of the observed Kui per Belt (|Bernstein et al.l 
120041 : iFuentes fc HolmanI 120081 ) and to date there are no 
detected Kuiper belt objects with these inclinations. In- 
tuitively it seems unlikely that we would have failed to 
detect this massive of a population, implying that the 
primordial gas disk must have had a small radial extent 
to avoid producing this population in the first place. 

To emphasize the differences between our hypothetical 
populations and the detected population we show the 
perihelion distance versus the inclination for all stable 
particles produced in all simulations as compared to the 
detected KBOs in figure [6l The observed Kuiper Belt 
occupies a distinctly different region of parameter space 
than any of the populations produced here. We note that 
there are a few objects detected on the periphery of the 
area of parameter space that is populated by these mod- 
els; for instance 2005 NU125 (i = 56.5°, a = 44.2 AU, 
e = 0.026). However, fewer than 1% of particles pro- 
duced in simulations with 100 AU disks have inclinations 
as low as 56.5°, and fewer than 5% of particles produced 
by 80 AU disks have inclinations this low. Therefore, 
if this object is a member of the population described 
in this paper we would expect to have detected 20-100 
KBOs with higher inclinations. Additionally the smaller 
50 AU disks produce a few particles at lower inclination 
and smaller perihelia that could potentially be consistent 
with observed KBOs. However these are still relatively 
low probability particles (< 10%) and, as we discuss later 
in this section, the current data sets place limits on these 
50 AU size disks. 

Using t wo surveys, the Deep Ecliptic Survey 
(DES; .Millis et al.l l2002[ ) and the Palomar Survey 



(jSchwamb et al.l 120101 ). we can quantify this statement 
and make rigorous observational limits. The DES sur- 
vey observed ov er 800 deg^ with a mean limiting VR 
magnitude 23.6 (lElliot et al.ll2005f ). Although this sur- 
vey targeted the ecliptic (±6°), it was still sensitive to 
these high inclination objects during a fraction of their 
orbits. The Palomar survey is a shallower survey (mean 
limiting VR magnitude 21.3) but with a wider coverage 
of 12,000 deg2 up to ±40° ecliptic latitude. Using the 
a,e, and i from our population of stable particles, we 
generate Attest = 10^ to 10^ test particles and create a 
synthetic population of orbits by assigning randomized 
orbital angles. For each of these test particles we deter- 
mine the probability (pij) that each particle (indexed by 
i) will be detected on each survey field (indexed by j) as 
a function of its intrinsic brightness {H). If a particle 
does not land on the field then pij{H) = 0. If it does 
land on the field, we can define a critical absolute mag- 
nitude (i^crit), a function of both the orbital properties 
and detection limit on the observation. 



crit,i,j 



^ miimitj - 2.5 log 



QiXi) 



(7) 



In this expression di is the distance between the sun and 
object, and is the distance between the observer and 
object (both measured in AU) and miimit,j is the back- 
ground limiting magnitude for the field. We parameterize 
the phase functi on (q(\i)) of these small bodies with the 
H G formalism (jBowell et al.lll989f ). ie. 



( rx\^-^^ 

^i(x)=exp ( -3.33 tan (^-j 



(8) 



02 (x) =exp ( —1.87 tan 



(f) 



1.22 



where Xi is the phase angle and G is a parameter 
describing the strength_of the o p position spik e (e.^ . 
Gehrels 1956: H aDke lHill [1991 iNelson et al.l [l998l : 
iRabinowitz et all [2007) However due to the fact that 
the fields are always chosen to view objects near opposi- 
tion, even varying this parameter from G=0 to 1 changes 
the critical value by less than 0.2 magnitudes so it does 
not notably affect our conclusions. While the probabil- 
ity of detection is a smoothed step function around this 
critical magnitude (i^crit), due to the uncertainties in 
the population size distribution, a simple step function 
will be adequate for these estimates. Therefore we use 
Pij{H) = 1 if H < i^crit,i,j, and zero otherwise. Addi- 
tionally, in order to determine the orbital elements of the 
detected object it must be detected in at least two dif- 
ferent observations. Therefore we define a new variable 
^(j?/) which we set to 1 if the two observations are tem- 
porally separated such that that the observed motion of 
the object is large enough to be resolvable by the instru- 
ment and small enough that the observations are clearly 
of the same object, and = otherwise. The total 

number objects we would expect to detect is then 



1 



Ntest A/'fields A/'fields pH 



detect " 



x5{j,j'){l-Usc,i)n{H)dH (9) 
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where fob 

i is the smah fraction of the field lost due 
to stellar obscuration, n{H) = dN/dH is the luminosity 
function for these hypothetical objects, and i^min and 
^max are the intrinsic brightness range for these parti- 
cles. 

To continue we must make a number of assump- 
tions about the properties of these hypothetical high- 
inclination objects. For this analysis we assume that the 
initial p lanetesimals have a K uiper belt-like size distri- 
bution. iBernstein et al.l (|2004|) found that the luminos- 
ity function can be well described by a combination of 
power-laws 

n{H) = no(10-^^ + clO-^'^)-\ (10) 

where /3 is the slope on the bright end, is the slope 
on the faint end, and c is chosen so that the break be- 
tween these two populations matches observations. Ob- 
servations indicate at least two distinctly different pop- 
ulations in the outer solar system. The dynamically hot 
Kuiper belt objects (KBOs) have a fairly shallow slope 
to their bright-end distribution, while the dynamically 
cold KBOs and the Jupiter Trojans have a steeper size 
distribution. The different origins of these size distribu- 
tions are not well understood. There is evidence that 
the hot population may have accumulated closer to the 
sun ^Levison fc SternI I2QQ1I : iGomesI I2QQ3I : iLevison et all 
l2QQ8f ). where the protoplanetary disk had shorter accre- 
tion times. For this reason we strongly prefer the as- 
sumption that the luminosity distribution of our hypo- 
thetical population is similar to that of the hot popula- 
tion, but we consider both for completeness. We revisit 
this line of reasoning at the end of this sect ion. 

Wh en looking at the bright end slope, iFraser et al.l 
(|2QlQf ) found different values for this power-law index /3 
in different populations, /3coid = 0.8, and /3hot = 0.35. 
Howe ver, the data cannot constrain the shallow end 
slope. IBernstein et al.l (|2004l ) found a break in both the 
hot and cold populations. We include this break in the 
cold populat i on and use /^coid = 0-38. However, we follow 
iFraser et al.l ()2010f ) and do not change the slope of the 
hot population after the "break" so /3{^^^ = 0.35. 

We assume that particles extend from 30 km to 1000 
km in size with a break at I) = 100 km. Additionally, 
we assume that the particle size distribution extends to 
1 km for the estimates of the mass, but we note that 
the calculation is insensitive to this choice as most of the 
mass lies in larger particles. To convert the number of 
expected detection into an estimate of the mass in this 
population we assume a bond albedo p = 0.05 and a 
density of p = 2 gcm~^. 

Finally, we can create a synthetic population and use 
the previously described method to calculate the ex- 
pected number of detections. For example, in our fiducial 
disk model, Rd = 100 AU, 7 = —1.5, and B = Ocut the 
fractional trapping efficiency is greater than 2 x 10~^. If 
the scattered particles have the same size distribution of 
the hot population, we expect the two surveys to have 
been able to detect ~ 10~^ of the scattered objects. Now, 
if we assume that the initial population of planetesimals 
was 40 M0, then, using the same size distribution, there 
would have been initially 10^ planetesimals greater that 
D > 30 km. Therefore if our solar system had had a 100 
AU disk gaseous disk we would expect the these two sur- 



veys to have detected over 100 objects. However, this 
number depends upon the initial assumption that we 
start out with a population of 40 planetesimals, a 
conservative lower limit. 

So perhaps a better way to think about this is to turn 
the argument around and ask, for the different disk sizes, 
what is the upper limit on this initial mass. Using a bino- 
mial distribution and a non-informative, Bayes-Laplace 
uniform prior we can determine the limit such that we 
are 95% confident that we will detect at least one par- 
ticle. For the populations used here, that corresponds 
to A^detect ^ 3. The results are plotted in figure [71 The 
upwards triangles indicate the results if we assume a size 
distribution consistent with the hot population (our pre- 
ferred model), while the open downwards triangles indi- 
cate the limits assuming the size distribution of the cold 
population. The error bars here are the statistical error 
bars (propagated through from figure [4|) . For disks that 
do not produce stable high-inclination particles we can 
not place any limits, and these disks are indicated with 
upwards arrows. 

If we assume that the size distribution of scattered 
planetesimals was similar to that of the hot population 
then we can rule out 80 and 100 AU disks as the limit on 
the initial population mass is below the expected theo- 
retical amount. Indeed, for disks with shallower profiles 
we can marginally rule out 60 AU disks with this size dis- 
tribution. However, if we relax this assumption on the 
size distribution then we cannot place as stringent limits 
on the disks with steeper surface density profiles. 

In order to significantly improve these constraints we 
would need a survey with a limiting magnitude similar to 
that of the DES (r=23.5) but with an area comparable 
to the Palomar survey (10,000 deg^). A non-detection in 
this type of survey would allow us to rule out all disks 
larger than 60 AU, except the steep exponential profile. 
The LSST, designed to cover 10, 000 deg^ with a typi- 
cal limiting magnitude of r=24.5 (jivezic et al.ll2008f ) will 
easily exceed this limit, allowing us to make meaningful 
constraints down to the 50 AU level. 

4.1. Discussion of Size Distribution 

In this paper we argue that the non-detection of this 
hypothetical high-inclination population can be used to 
constrain the size of the protosolar nebula. Of course, 
another possibility is that the population does exist, but 
that it has such a steep size-distribution that most of 
the mass is in bodies well below 100 km in size, and 
thus undetectable. There are three main reasons why we 
believe this is unlikely. 

First, there is another potential reservoir in which some 
remnants from this same initial population of planetes- 
imals are beheve to have been stored: the Sedna popu- 
lation. The 1,000 km Sedna is so far unique in that it 
has a perihelion distance of 76 AU, well detached from 
Neptune, but it has a semimajor a xis around 500 AU , 
well interior to the main Oort cloud (jBrown et al.ll200"3 ). 
The observational difficulties in detecting an object on 
this type of orbit implies there must be on the order of 
100 similar sized planetoids on similar orbits. The best 
dynamical mechanism to produce this type of orbit is 
to have an object scattered by a planet very early on 
in the solar system (at the same time considered in this 
paper) so that its pericenter can be lifted via an interac- 
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Fig. 7. — Observational constraints on the upper limit on the ini- 
tial mass of planetesimals in the planet forming region. As in figure 
[31 the left panel is for sharp cutoff disks and the right panel is for 
disks with exponential profiles. The color indicates the disk surface 
density profile (On the left, 7 = 1,3/2 in blue and black, respec- 
tively. On the right, 7 = 1/2,1 in blue and black, respectively.) 
Upwards triangles indicate results for a shallow size distribution 
consistent with the hot population, while the downwards triangles 
indicate the steeper profile consistent with the cold population (see 
text). Arrows indicate disk parameters for which we cannot limit 
the initial mass. 

tion with the birth cluster potential or by the passage of 
nearby stars (iMorbidelli & LevisonI |2QQ5 iBrasser et alJ 
[2OO6). The large predicted number of 1,000 km-sized ob- 
jects implies a shallow size distribution, and constraints 
from the TAOS survey indicate th at this population must 
have P < 0.88 (jWang et alJl2009[ ). i.e. the size distribu- 
tion cannot be much steeper than that of the cold popu- 
lation. 

Second, in the Kuiper belt where we can directly mea- 
sure the size distributions, we see that the hot population 
has a shallower size distribution than the cold popula- 
tion. As mentioned previously, there is evidence that 
the hot population may have accumulated closer to the 
sun jLevison fc Sternll2QQl IGomesl[2QQl iLevison et all 
|2QQ8[ ). where the protoplanetary disk had shorter accre- 
tion times. Therefore we would expect that the source 
region for our hypothetical population, being even closer 
to the sun than the hot population, may have a shallower 
population. At least it seems unlikely to be significantly 
steeper. 

Finally, we can be confident that the population that 
produced the comets could not have had a very steep size 
distribution. If it had, then mutual collisions would have 
ground the comets down bef ore they could have bee n 
scattered into the Port Cloud (|Stern fc Weissma"nll2QQl[ ). 
iCharnoz fc Morbidellil (|2QQ3f ) demonstrated that this col- 
lisional grinding problem can be avoided by assuming a 
shallower size distribution up through 10s to 100s of km, 
such as the size distributions used in this paper. While 
the Oort cloud likely has a source region exterior to the 
population considered here (the Nice disk,lTsigan is et al.l 
|2005) or even from other stars (|Levison et al. 2010). the 
same argument holds for the Sedna region. For these 
reasons we believe that the size distribution of this pop- 
ulation are unlikely to have a steeper distribution than 
the cold population and very likely have been similar to, 
or even more shallow than, the hot population. 

5. CONCLUSIONS AND DISCUSSION 

In this paper we present a dynamical mechanism by 
which high inclination Kuiper Belt objects can be pro- 
duced during the planet formation process due to Kozai 



cycles excited by the presence of the gaseous protoplan- 
etary disk. If our own solar system's disk had been ra- 
dially extended {Rd ^ 50 AU) a significant number of 
scattered planetesimals would have been placed onto high 
inclination orbits which would be stable for the lifetime 
of the solar system. We find that the current observed 
lack of objects on these orbits in the DES and Palomar 
Distance Solar System Survey suggests that we can ex- 
clude very extended disks 80 AU) so long as we are 
correct in our assumption of a size distribution more sim- 
ilar to the hot KBOs than the cold KBOs. We argue that 
while we cannot technically rule out other options such 
as that the planetesimals in the planet forming region 
were orders of magnitude less numerous than expected 
by theory or that the size distribution was dominated by 
small planetesimals, we believe these options to be quite 
unlikely. 

Future observations from deep all sky surveys (such 
as the LSST) will allow us to place more stringent con- 
straints on the size of the disk. With this type of sur- 
vey, if we are able to detect a population of Kuiper belt 
objects with i ~ 60° the orbital distribution of these 
particles will reveal information on the shape of the pro- 
toplanetary disks, especially revealing information about 
the truncation. However, it is also possible (and perhaps 
likely) that future observations will not reveal this pop- 
ulation. If this is the case we will be able to conclude 
that the disk must have been smaller than 50 AU, but 
we will not be able to place more stringent limits because 
if the disk had been smaller than 50 AU, the dynamical 
mechanism presented in this paper cannot produce stable 
particles. 

This method if the first technique to constrain the size 
of the gaseous protoplanetary disk of a known plane- 
tary system during the time of planet formation. The 
disk sizes we find consistent with observations implies 
that in some ways our solar system is typical. Although 
many of the first protoplanetary disks observed were 
quite radially extended, smaller disks with sizes below 60 
AU are quite commo n according to recent observations 
([Andrews et al.|[2QTQl ). Additionally, this small disk size 
is consistent with other pictures of the early solar system. 
This work is consistent with the idea that a sharp cut- 
off in the planetesimal disk at around 50 AU might have 
corresponded to a sharp cut off in the gas disk. Future 
deep all-sky surveys will allow us to probe closer to this 
50 AU radius and allow us to meaningfully constrain on 
the theoretical models of of solid-gas interactions during 
the time of the formation of the solar system. 



We thank M. Schwamb for guidance in our comparison 
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